* 
* Replication file for 
* "Leave It As It Is": International Network Effects on Protected Lands
* Darren Hawkins and Jay Goodliffe
* International Interactions
*

* install required extension:
ssc install reghdfe

* set working directory 
cd "directory path here"

* load dataset
use "protectedlands.dta", clear

* set dataset for TSCS analyses
xtset country year, yearly


* General Approach: Lagged Dependent Variable and Country Fixed Effects with PCSEs

* Table 1

* domestic variables only 
xtpcse perc_protected_bio polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
* country fixed effects jointly significant?
testparm i.country
* domestic variables jointly significant?
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
* xtpcse does not estimate log-likelihood, so does not calculate BIC
* use areg (or xtreg), which produces same coefficients and log-likelihood (but different SEs)
quietly areg perc_protected_bio polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic

* trade network only
xtpcse perc_protected_bio L.tradeprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
quietly areg perc_protected_bio L.tradeprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic

* IO network only
xtpcse perc_protected_bio L.igo10protected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
* joint significance of subset of domestic variables
test polity2 lngdppc perc_rural_pop
quietly areg perc_protected_bio L.igo10protected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic

* aid network only
xtpcse perc_protected_bio L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
test polity2 lngdppc perc_rural_pop
quietly areg perc_protected_bio L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic

* regional diffusion only
xtpcse perc_protected_bio L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
test polity2 lngdppc perc_rural_pop
quietly areg perc_protected_bio L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic

* full model: all networks and region
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
* joint significance of IO and aid networks
test L.igo10protected L.aidprotected
test polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
test polity2 lngdppc perc_rural_pop
quietly areg perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic
* save this sample for the next model
capture drop sample
gen sample = e(sample)

* reduced model: drop region, density and ruleoflaw and natrents
* has best BIC
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc perc_rural_pop L.perc_protected_bio i.country if sample, pairwise
testparm i.country
test L.igo10protected L.aidprotected
test polity2 lngdppc perc_rural_pop
quietly areg perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc perc_rural_pop L.perc_protected_bio if sample, absorb(country)
estat ic

* Fn 6 (Also SI p 12): other operationalizations of IO network
* linear power transformation (ie, no transformation)
xtpcse perc_protected_bio L.tradeprotected L.igoprotected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
* quadratic power transformation
xtpcse perc_protected_bio L.tradeprotected L.igosqprotected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
* importance-weighted 
xtpcse perc_protected_bio L.tradeprotected L.igowtprotected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
* GDP-weighted
xtpcse perc_protected_bio L.tradeprotected L.igogdpwdiusprotected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise


* Figure 10: substantive effects
* 2 sd * long-range effect [2 sd * b(var) / (1-b(LDV))]
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
capture drop sample
gen sample = e(sample)
xtsum perc_protected_bio if sample
xtsum L.tradeprotected if sample
scalar trade_sd = r(sd_w)
xtsum L.igo10protected if sample
scalar io_sd = r(sd_w)
xtsum L.aidprotected if sample
scalar aid_sd = r(sd_w)
xtsum L.noaid if sample
scalar noaid_sd = r(sd_w)
xtsum L.regionprotected if sample
scalar region_sd = r(sd_w)
xtsum polity2 if sample
scalar polity_sd = r(sd_w)
xtsum lngdppc if sample
scalar gdppc_sd = r(sd_w)
scalar gdppc_m = r(mean)
xtsum lnnatrents if sample
scalar natrents_sd = r(sd_w)
xtsum lnpopdensity if sample
scalar popdensity_sd = r(sd_w)
xtsum perc_rural_pop if sample
scalar rural_sd = r(sd_w)
xtsum rle if sample
scalar rle_sd = r(sd_w)
nlcom (IO: 2*io_sd*_b[L.igo10protected]/(1-_b[L.perc_protected_bio])) ///
	(noaid: 2*noaid_sd*_b[L.noaid]/(1-_b[L.perc_protected_bio])) ///
	(trade: 2*trade_sd*_b[L.tradeprotected]/(1-_b[L.perc_protected_bio])) ///
	(aid: 2*aid_sd*_b[L.aidprotected]/(1-_b[L.perc_protected_bio])) ///
	(region: 2*region_sd*_b[L.regionprotected]/(1-_b[L.perc_protected_bio])) ///	
	(polity: 2*polity_sd*_b[polity2]/(1-_b[L.perc_protected_bio])) ///
	(natrents: 2*natrents_sd*_b[lnnatrents]/(1-_b[L.perc_protected_bio])) ///
	(popdensity: 2*popdensity_sd*_b[lnpopdensity]/(1-_b[L.perc_protected_bio])) ///
	(rle: 2*rle_sd*_b[rle]/(1-_b[L.perc_protected_bio])) ///
	(rural: 2*rural_sd*_b[perc_rural_pop]/(1-_b[L.perc_protected_bio])) ///
	(gdppc: 2*gdppc_sd*_b[lngdppc]/(1-_b[L.perc_protected_bio])), post
parmest, format(estimate min95 max95 %8.4f) saving("substeffects.dta", replace)

use "substeffects.dta", clear
egen varnum = fill(12/8 6/1)
twoway (rspike min95 max95 varnum, horizontal) (scatter varnum estimate, mcolor(dknavy)), ///
	ylabel(1 "GDP/capita" 2 "Rural Population" 3 "Rule of Law" 4 "Population Density" 5 "Resource Rents" ///
	6 "Democracy" 7 " " 8 "Regional Diffusion" 9 "Aid Network" 10 "Trade Network" 11 "Receives No Aid" 12 "IO Network", angle(horizontal) labsize(medium)) ///
	title("") ytitle("") xline(0, lcolor(dknavy) lpattern(dot))  ///
	legend(off) ytick(1/6 8/12) xlabel(, format(%8.2f)) graphregion(color(white))
graph save "substeffects", replace
graph export "substeffects.pdf", replace


* Supplemental Information

* Table SI-1: Summary Stats
use "protectedlands.dta", clear
* establish dataset used
areg perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
capture drop sample
gen sample = e(sample)
summarize perc_protected_bio L.tradeprotected L.igo10protected L.noaid L.regionprotected polity2 ///
	gdppc_2010 lngdppc natrents lnnatrents popdensity lnpopdensity perc_rural_pop rle if sample, d
* aid network measured only 
summarize L.aidprotected if noaid == 0 & sample == 1, d
* get the within sd
xtsum perc_protected_bio L.tradeprotected L.igo10protected L.noaid L.regionprotected polity2 ///
	gdppc_2010 lngdppc natrents lnnatrents popdensity lnpopdensity perc_rural_pop rle if sample
xtsum L.aidprotected if noaid == 0 & sample == 1


* Table SI-2: Alternative temporal dependence
use "protectedlands.dta", clear
* original full model
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test L.igo10protected L.aidprotected
* Prais-Winsten (AR1 instead of LDV)
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle i.country, pairwise c(ar1)
testparm i.country
* country-specific AR1 (instead of LDV)
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle i.country, pairwise c(psar1)
testparm i.country
matrix rhos = e(rhomat)'
matlist rhos
* Cut and paste into Excel, and find median
* LDV AND AR1 (Prais-Winsten)
* trade, IO, noaid (aid); polity
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise c(ar1)
testparm i.country


* SI fn 4: serial correlation/unit root
* Easiest to do with Correia (2017) estimator installed above (reghdfe)
* unit root
reghdfe perc_protected_bio L.perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country) resid(ehatldvfe)
xtunitroot fisher ehatldvfe, dfuller lags(0)
* reject unit root in the residuals
* any serial correlation in 
regress ehatldvfe L.ehatldvfe L.perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
display e(N)*e(r2)
display 1-chi2(1, e(N)*e(r2))
* statistically significant p = 0.026, but substantively small (rho = 0.033)
* try 2 lags on DV
reghdfe perc_protected_bio L(1/2).perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country) resid(ehatldv2fe)
regress ehatldv2fe L.ehatldv2fe L(1/2).perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle
display e(N)*e(r2)
display 1-chi2(1, e(N)*e(r2))
* made it worse: rho = -0.081


* SI p 6: compare different lags using BIC
* SEs do not affect IC, so just use reghdfe
* 1 > 0&1 (ie, BIC1 < BIC0&1)
reghdfe perc_protected_bio L.perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country) 
estat ic
reghdfe perc_protected_bio L(0/1).tradeprotected L(0/1).igo10protected L(0/1).aidprotected L(0/1).noaid L(0/1).regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
estat ic
* 1 > 1&2
* need to make sure run on same dataset/observations
reghdfe perc_protected_bio L.perc_protected_bio L(1/2).tradeprotected L(1/2).igo10protected L(1/2).aidprotected L(1/2).noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country)
gen sample2 = e(sample)
estat ic
reghdfe perc_protected_bio L.perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle if sample2, absorb(country) 
estat ic
drop sample2
* 1 > 2
reghdfe perc_protected_bio L2.tradeprotected L2.igo10protected L2.aidprotected L2.noaid L2.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
gen sample2 = e(sample)
estat ic
reghdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio if sample2, absorb(country)
estat ic
drop sample2
* 1 > 3
reghdfe perc_protected_bio L3.tradeprotected L3.igo10protected L3.aidprotected L3.noaid L3.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country)
gen sample3 = e(sample)
estat ic
reghdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio if sample3, absorb(country)
estat ic
drop sample3

* SI fn 6: dynamics quantities of interest
* same for any variable
use "protectedlands.dta", clear
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
scalar error_corr_rate = 1-_b[L.perc_protected_bio]
scalar long_run_effect = ((_b[L.tradeprotected])/(1-_b[L.perc_protected_bio]))
scalar period1 = _b[L.tradeprotected]
scalar period1prop = period1/long_run_effect
scalar p1cumprop = period1prop
forvalues i = 2/20 {
	local im1 = `i' - 1
	scalar period`i' = _b[L.perc_protected_bio]*period`im1'
	scalar period`i'prop = period`i'/long_run_effect
	scalar p`i'cumprop = p`im1'cumprop + period`i'prop
	}
forvalues i = 1/20 {
	display "period `i' effect = " period`i'
	}
forvalues i = 1/20 {
	display "period `i' effect as a proportion of long run effect = " period`i'prop
	}
forvalues i = 1/20 {
	display "cumulative period `i' effect as a proportion of long run effect = " p`i'cumprop
	}
* median lag: 8 years
display "error correction rate = " error_corr_rate



* Table SI-3: SEs and FEs
use "protectedlands.dta", clear
* original full model
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
testparm i.country
test L.igo10protected L.aidprotected
* country-clustered SEs
* Correia (2017) estimator installed above
reghdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio, absorb(country) vce(cluster country)
test L.igo10protected L.aidprotected
* need to run as regular regression to test FEs
reg perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, vce(cluster country)
testparm i.country
* Year FEs instead of LDV (country-clustered SEs)
reghdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country year) vce(cluster country)
reg perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle i.country i.year, vce(cluster country)
testparm i.country
testparm i.year


* Table SI-4: Functional Form
use "protectedlands.dta", clear
* exponential model with country FE and LDV
* logging the DV changes from right-skewed to left-skewed
* create DV to be lagged
gen logprotect = ln(perc_protected_bio)
* fix the zeros close to the lowest non-zero
replace logprotect = -14 if perc_protected_bio == 0
* Correia SEs
ppmlhdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.logprotect, absorb(country) vce(cluster country)
test L.igo10protected L.aidprotected
* p = 0.06
* To test joint significance, use poisson w/FEs
poisson perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.logprotect i.country, vce(cluster country)
testparm i.country

* panel fractional probit
* Papke and Wooldridge estimator uses fracreg, country-variable averages instead of FEs, and year FEs instead of LDV
* establish sample
capture drop sample
reghdfe perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle, absorb(country) vce(cluster country)
gen sample = e(sample)
* create country-variable averages
* cannot use the L. approach with egen, by
gen trade_1 = L.tradeprotected
gen IO_1 = L.igo10protected
gen region_1 = L.regionprotected
gen aid_1 = L.aidprotected
gen noaid_1 = L.noaid
egen mn_trade = mean(trade_1) if sample, by(country)
egen mn_IO = mean(IO_1) if sample, by(country)
egen mn_aid = mean(aid_1) if sample, by(country)
egen mn_noaid = mean(noaid_1) if sample, by(country)
egen mn_region = mean(region_1) if sample, by(country)
egen mn_polity = mean(polity2) if sample, by(country)
egen mn_gdppc = mean(lngdppc) if sample, by(country)
egen mn_natrents = mean(lnnatrents) if sample, by(country)
egen mn_popdens = mean(lnpopdensity) if sample, by(country)
egen mn_rural = mean(perc_rural_pop) if sample, by(country)
egen mn_rle = mean(rle) if sample, by(country)
fracreg probit perc_protected_bio trade_1 mn_trade IO_1 mn_IO aid_1 mn_aid noaid_1 mn_noaid region_1 mn_region polity2 mn_polity lngdppc mn_gdppc lnnatrents mn_natrents lnpopdensity mn_popdens perc_rural_pop mn_rural rle mn_rle i.year, vce(cluster country)
* joint significance tests
testparm mn_*
testparm i.year
* now get APEs with bootstrapped clustered SEs
preserve
keep if sample
xtset, clear
capture program drop APE_boot
program APE_boot, rclass
fracreg probit perc_protected_bio trade_1 mn_trade IO_1 mn_IO aid_1 mn_aid noaid_1 mn_noaid region_1 mn_region polity2 mn_polity lngdppc mn_gdppc lnnatrents mn_natrents lnpopdensity mn_popdens perc_rural_pop mn_rural rle mn_rle i.year, vce(cluster country)
predict x1b1hat, xb
gen scale=normalden(x1b1hat)
gen petrade=scale*_b[trade_1]
summarize petrade
return scalar apetrade=r(mean)
gen peIO=scale*_b[IO_1]
summarize peIO
return scalar apeIO=r(mean)
gen peaid=scale*_b[aid_1]
summarize peaid
return scalar apeaid=r(mean)
gen penoaid=scale*_b[noaid_1]
summarize penoaid
return scalar apenoaid=r(mean)
gen peregion=scale*_b[region_1]
summarize peregion
return scalar aperegion=r(mean)
gen pepolity=scale*_b[polity2]
summarize pepolity
return scalar apepolity=r(mean)
gen pegdppc=scale*_b[lngdppc]
summarize pegdppc
return scalar apegdppc=r(mean)
gen penatrent=scale*_b[lnnatrents]
summarize penatrent
return scalar apenatrent=r(mean)
gen pepopden=scale*_b[lnpopdensity]
summarize pepopden
return scalar apepopden=r(mean)
gen perural=scale*_b[perc_rural_pop]
summarize perural
return scalar aperural=r(mean)
gen perle=scale*_b[rle]
summarize perle
return scalar aperle=r(mean)
drop x1b1hat scale petrade peIO peaid penoaid peregion pepolity pegdppc penatrent pepopden perural perle
end
* APEs and Bootstrapped SE within countries
bootstrap r(apetrade) r(apeIO) r(apeaid) r(apenoaid) r(aperegion) r(apepolity) r(apegdppc) r(apenatrent) r(apepopden) r(aperural) r(aperle), reps(10000) seed(123456) cluster(country) idcluster(newid): APE_boot
program drop APE_boot
restore

* Stock & Watson approach to exogeneity using Newey-West 
* the L?D shows the cumulative effect up to ? and the L? shows the cumulative effect at ?
* tried this with lags up to 8: long-term and short-term effects simliar
* for other lags of length n, replace 2 with n-1, and 3 with n in the foreach loop
foreach var of varlist tradeprotected igo10protected aidprotected noaid regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle {
	forval l = 1/2 {
		gen L`l'D_`var' = L`l'.D.`var'
	}
	gen L3_`var' = L3.`var'
}
newey perc_protected_bio L?D_tradeprotected L?_tradeprotected L?D_igo10protected L?_igo10protected L?D_aidprotected L?_aidprotected L?D_noaid L?_noaid L?D_regionprotected L?_regionprotected L?D_polity2 L?_polity2 L?D_lngdppc L?_lngdppc L?D_lnnatrents L?_lnnatrents L?D_lnpopdensity L?_lnpopdensity L?D_perc_rural_pop L?_perc_rural_pop L?D_rle L?_rle i.country if sample, lag(3) force
* count number of countries included
egen numc = group(country) if e(sample)
sum numc
drop numc

* SI fn 11: longer se lag
newey perc_protected_bio L?D_tradeprotected L?_tradeprotected L?D_igo10protected L?_igo10protected L?D_aidprotected L?_aidprotected L?D_noaid L?_noaid L?D_regionprotected L?_regionprotected L?D_polity2 L?_polity2 L?D_lngdppc L?_lngdppc L?D_lnnatrents L?_lnnatrents L?D_lnpopdensity L?_lnpopdensity L?D_perc_rural_pop L?_perc_rural_pop L?D_rle L?_rle i.country if sample, lag(12) force


* SI p 12: developed vs developing
use "protectedlands.dta", clear
* developing
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country if developed == 0, pairwise
* developed
xtpcse perc_protected_bio L.tradeprotected L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country if developed == 1, pairwise


* SI p 12: dependence on trade 
* no difference in those with more trade exposure
xtpcse perc_protected_bio c.L.tradeprotected##c.tradepropgdp L.igo10protected L.aidprotected L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
test tradepropgdp L.tradeprotected#tradepropgdp


* SI p 12: dependence on aid
* no effect of absolute level of aid (interaction with aid network)
xtpcse perc_protected_bio L.tradeprotected L.igo10protected cL.aidprotected##c.aid_millions L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise
* no effect of Aid/GDP (interaction with aid network)
xtpcse perc_protected_bio L.tradeprotected L.igo10protected cL.aidprotected##c.aidpropgdp L.noaid L.regionprotected polity2 lngdppc lnnatrents lnpopdensity perc_rural_pop rle L.perc_protected_bio i.country, pairwise

